Seasonal variation in wing size and shape of Drosophila melanogaster reveals rapid adaptation to environmental changes

Populations in seasonal fluctuating environments receive multiple environmental cues and must deal with this heterogenic environment to survive and reproduce. An enlarged literature shows that this situation can be resolved through rapid adaptation in Drosophila melanogaster populations. Long-term monitoring of a population in its natural habitat and quantitative measurement of its responses to seasonal environmental changes are important for understanding the adaptive response of D. melanogaster to temporal variable selection. Here, we use inbred lines of a D. melanogaster population collected at monthly intervals between May to October over a temporal scale spanning three consecutive years to understand the variation in wing size and wing shape over these timepoints. The wing size and shape of this population changed significantly between months and a seasonal cycle of this traits is repeated for three years. Our results suggest that the effects of environmental variables that generated variation in body size between populations such as latitudinal clines, are a selective pressure in a different manner in terms of seasonal variation. Temperature related variable have a significant nonlinear relation to this fluctuating pattern in size and shape, whereas precipitation and humidity have a sex-specific effect which is more significant in males.


Material and methods
The wild population of Drosophila melanogaster was collected from an orchard in Yeşilöz, Turkey (40.30° N, 32.34° E) monthly from June to October in 2014 and May to October in 2015 and 2016. Every collection occurred approximately with a 30-day interval with traps hung always in the same location for 1-2 days. Isofemale lines were established upon collection and were maintained at standard conditions; 25 °C and 60-65% humidity on a 12:12 h light/dark cycle on a standard cornmeal-agar-sugar-yeast medium.
After two generations, inbred lines were established from isofemale lines through 20 generations of full-sib mating which is considered in a theoretically inbreeding level of f ≅ 1 70 . A minimum of 10 inbred lines per collection month and year were used for the wing measurements, summing up to a total of 166 inbred lines in this study. One exception to this is June 2014 sampling time-point, which has been phenotype only for 6 lines due to low survival upon inbreeding.
Wing measurements. To obtain flies for the wing size and shape measurements, we put 25 adult flies of each sex in 3 replicate bottles containing 50 mL of standard medium. To prevent overcrowding the parents were discarded after 3 days from the bottles. Upon eclosion, at least 25 flies of each sex moved into egg-collection chambers with %2 agar medium, and a pile of yeast paste at 25 °C. After 6-8 h, 150 eggs are collected from the surface with fine forceps and transferred evenly to 3 replicate vials containing 7 mL standard medium. The collected eggs are allowed to develop under standard conditions, as described above. Three days after eclosion adult flies from every replicate vial were collected and stored in ethanol. These standard conditions are used for minimizing the environmental effects on body size. www.nature.com/scientificreports/ A minimum of 15 individuals per sex have been dissected from 166 inbred lines, summing up to a total of 5377 wings being measured. The left wing was removed with fine forceps and mounted on a microscope slide in Entellan® (Merck Millipore). All wings were photographed by using a camera attached stereo microscope (Leica S9i) and then digitized as jpeg format.
Analysis. For the estimation of wing shape, 11 landmarks were digitized of the left wing ( Fig. 1) of each fly using tpsUtil and tpsDig2 71,72 . The tps files were imported to the MorphoJ software version 1.07a 73 . To remove the position and orientation variations from the coordinates the raw landmark coordinates were aligned and superimposed using a Procrustes Fit function and wing shape variation was analyzed. The effects of sex and sampling months on shape were tested with Procrustes ANOVA. MANOVA was performed on shape coordinates (Procrustes coordinates) to test whether the population diverged in shape by months. To maximize the differences between groups by comparing mean shapes a Canonical variate analysis (CVA) 73 was used to analyze discrimination between groups where month used as a classifier. The pairwise differences in shape were analyzed by using a permutation test (10.000 rounds) with Procrustes distances, which is the sum of the distance between the landmarks. Additionally, a discriminant function analysis (DFA) was used to find the shape differences between sampling months. To test for allometry (relationship of body size to shape) a multivariate regression of shape (the Procrustes coordinates as dependent variable) on size (centroid size as independent variable) was performed, pooled within sub-groups of lines, with a permutation test with 10,000 rounds for testing statistical significance. All analyses were performed using MorphoJ version 1.07a 73 .
Wing centroid sizes (computed from raw data of landmarks) were used as a measure of wing size. Male and female flies were analyzed separately due to the significant sexual dimorphic differences in body size. The centroid size (CS) data were first analyzed by a Shapiro-Wilk normality test. As the data did not follow a normal distribution, analyses were continued by nonparametric tests. To determine the difference between the three replicate of wing measurements, a Kruskal-Wallis test was performed. Since there was no significant difference between replicas (p > 0.05), the data was pooled. Differences among months of CS were computed with the Kruskal-Wallis test followed by Dunn test with Bonferroni correction for multiple comparisons and pairwise tests using the "ggstatsplot" package 74 . All analyses were performed in R software 75 To test the relationship between climatic variables and centroid size we used climatic variables recorded monthly by the Turkish State Meteorological Service in Güdül, 6 km far from the sampling point (Table 1). Monthly mean temperature (T mean ), monthly mean maximum temperature (T max ), monthly mean minimum temperature (T min ), monthly mean humidity (H mean ), monthly mean maximum humidity (H max ), monthly mean minimum humidity (H min ) and total monthly precipitation (P total ) variables were used. The analyses were split in two parts: (i) using the climatic data for 15 days before the collection dates, (ii) using the climatic data for the collection days. Our preliminary results supported a strong relationship with climatic variables prior the collection dates and analyses were performed with these variables (given in Table 1).
We performed a nonparametric regression analysis using Generalized Additive Models (GAMs) 76 to estimate temporal trends in wing CS and its relationships to climatic variables. GAMs provide a capture of linear and non-linear relationships between variables 77 . Thus, the model allows better prediction of trends in the means of the variables and associated uncertainties. As the climatic categories were collinear, we have adjusted several models with a single independent variable. The dependent variable was CS and independent variables were monthly mean temperature (T mean ), monthly mean maximum temperature (T max ), monthly mean minimum temperature (T min ), monthly mean humidity (H mean ), monthly mean maximum humidity (H max ), monthly mean minimum humidity (H min ) and total monthly precipitation (P total ). Models were run using "mgcv" package (version 3.6.2) 78 in R 75 with residual marginal likelihood (REML) smoothness selection 78,79 . The corresponding R function is gam(CS ~ s(temperature related variable), method = 'REML', data) and gam(CS ~ s(humidity and precipitation related variable), method = 'REML', data). For CS analyses, final model selection was determined by calculating the Akaike Information Criterion (AIC) for all models where the type of smooth terms specify thin plate regression splines (default smooth for s terms) or cubic regression splines (bs = "cs"). www.nature.com/scientificreports/

Results
Wing size variation. We studied the effects of the monthly changing environments on wing size and shape in a temperate population from Yeşilöz, Turkey collected over successive years and months. Inbred lines from each collection month were used for size measurements. While size is a sexual dimorphic trait and CS differ significantly between sexes (χ 2 = 3985.8, p < 0.0001) (Fig. S1), male and females wing measurements were analyzed separately.
The mean and standard errors and coefficients of variation for wing measurements for each year and month were given in Table 2 and we found evidence of seasonally fluctuating wing size which is repeated over three years. Mean wing centroid size was smaller in May and rises from June to July, which is followed by smaller sizes in August to September then an increase in October ( Table 2, Fig. 2). This fluctuating pattern is repeating every year, except for October 2016, where the wing sizes are decreasing and that was not consistent with the pattern seen in previous years (Fig. 2). Inbred lines originating from June and July are one average bigger than other months in both sexes ( Table 2). The coefficients of variation were particularly higher in some months (Table 2). This might be due to a greater genetic heterogeneity between inbred lines in October 2014, June 2015, August 2015 and May 2016.  Fig. 2) as also in males (χ 2 = 144.5, df = 5, p-value < 0.0001, Fig. 2). Pairwise comparisons performed with Dunn test showed a significant difference between months except May-August, May-September, May-October, June-July and August-October in females. Whereas pairwise comparisons between months in males shows that only June and July differ significantly (p < 0.0001) from other months, all other comparisons are not significant (Fig. 2).
Wing shape variation. We digitized 11 landmarks to evaluate differences in wing shape between months.
The canonical variate analysis provided a clear discrimination of the months defined by the two first canonical axes which respectively accounted for 38.97% and 26.20% of the among-month variance in females (total 65.17%) and 35.40% and 25.36% of the among-month variance in males (total 60.76%). The results of the canonical variate analysis with Mahalanobis distance and Procrustes distances showed morphological variability across the months (Table S1). Permutation tests (10,000 rounds) for Procrustes distances and Mahalanobis distances among months were significant in all comparisons (Table S1) where May and October lines exhibited difference less than others. The highest Procrustes distance in shapes were found among May and July, where the lowest distance was found between May and October in both sexes. Furthermore, there was a significant difference between sexes (Mahalanobis distance among groups = 2.1648, p < 0.0001) within the population (Figs. S2, S3). The multivariate regression coefficient showed that the effect of size (centroid size) on shape (Procrustes coordinates) was highly significant (P < 0.0001). In addition, we observed a moderate amount of allometry that accounted for 22.06% of the total shape variance. In other words, line specific shape variation explained by size differences.
Climatic effects for size and shape of wing. To find the environmental factor that led to this repetitive pattern across years we applied a Generalized Additive Model (GAM) analysis to assess the relationship between wing centroid size and the climatic variables.
To estimate the relationship between the dependent and independent variables, the GAM analysis was used. The findings of these relationships were estimated as linear and nonlinear effects using smoothing splines (Fig. S4), results presented in Table 3. Results in Table 3 showed that some models resulted with the effective degrees of freedom (edf) values near to 1, the measures suggest a linear pattern. Whereas edf values were larger than 1 suggested a nonlinear pattern in the relationship between the dependent and independent variable. Furthermore, the model output has indicated that all spline terms of climate variables were significant variables of wing size (Table 3).
Temperature related variables have a nonlinear relationship with wing size. Temperature related variables fit well with a cubic regression spline (model selection made by AIC criterion). On the other hand, humidity has a linear effect (edf = 1) on male and females in terms of both mean and minimum values. However, precipitation has a nonlinear effect in females, but a linear effect in males (Table 3).

Discussion
Seasonal fluctuations are main sources for rapidly changing environments in temperate regions and recent studies provided evidence for rapid change in SNP frequencies among seasons [12][13][14] . Previous studies also revealed that seasonal variation in adaptive phenotypic traits like thermal and starvation tolerance 22 , and innate immunity 25 in D. melanogaster populations and all together suggesting that many loci respond to seasonal selection. This Table 3. Results of a GAM model, with the effective degrees of freedom (edf), F statistics, and P value using non-parametric smoothers. www.nature.com/scientificreports/ evidence raises new questions about seasonal adaptation and adaptive phenotypic traits. This is main reason to carry out this study to investigate the presence of seasonal variation in wing size and shape. We have investigated the seasonal variation in wing size and shape in a temperate Drosophila population. Our results show that wing size and shape changed significantly between months and the observed seasonal cycle is repeated for three consecutive years. We show that temporal selection acts along seasons in the Yeşilöz population of D. melanogaster and wing size and shape exhibit a rapid response among months. Our analysis reveals evidence for temporal varying selection on wing size and shape; and this suggesting differences within population between month from this geographic region. Body size is characterized by a complex genetic nature, that a large number of genes are involved, and the alleles have small effects of this trait. It is hard to find the adaptive polymorphism for such traits to towards more simple traits with alleles of large effects. Furthermore, the methodology used in this study provides measurements of inbred lines derived monthly from the population, reflecting a random genotypic sample in the population. Our results of the panel of inbred lines support that the source of seasonal variation in wing size and shape is an outcome of rapid seasonal adaptation, as the frequency of some alleles of different loci vary between months and contribute to the variation of size. It is likely that many ecological and environmental variables that changed between months are interacted to generate variation in body size and strong selection occurs between months for this morphological trait.
The wing size and shape variations depend on the amount of genetic variation in the population. The fluctuating response to selection shows that the genetic variation is maintained in the population but allele frequency for loci responsible for size and shape are changed rapidly through seasons. This might partially explain the difference within months and between inbred lines. For example, some individuals of an inbred line might have outliers within a month, however a general pattern persists. Previous studies showed that SNP frequencies in D. melanogaster populations vary seasonally [12][13][14] and transposable elements (TEs) were significantly associated with seasons 80 . Thus, it is possible that selection favors size and / or size related traits depended on changing environmental variables as our wing size and shape results support.
Temperature has also direct environmental effects on body size in Drosophila and other ectotherms and it is known as the major driver for body size clines 35,59,60,62,[81][82][83][84][85] . Our results indicated that wing size and shape significantly affected by monthly temperature, humidity and precipitation. Accelerated precipitation and humidity have a significant linear relationship with increased wing size especially on males (Table 3). In contrast to this, temperature has a significant non-linear relationship with wing size (Table 3). Our results show that seasonally varying precipitation and humidity have a higher selective pressure on body size than temperature. Increased humidity and precipitation promote larger body size significantly (Table 3, Fig. S4). Our results are not consistent with some other studies about the relationship between size and precipitation and humidity. For example, the wing size of Bactrocera tryoni, an Australian endemic horticultural pest species, is also affected by precipitation but shows a negative relationship with precipitation seasonality 86 . Stilwell et al. 40 proposed that the latitudinal cline in body size of the seed-feeding beetle Stator limbatus generated by humidity and seasonality then temperature, and the results support humidity for a better candidate to explain the body size clines where body size increased with decreased humidity. In a similar manner Gibb et al. 87 showed that the body size declines as precipitation increases in ant species. Experimental evolution at low and high humidity showed in larger wing areas at low humidity conditions in D. melanogaster 88 . For D. simulans and D. mercatorum, Przybylska et al. 89 found that the flies were larger in the dry season. Generally, the relationship between humidity and body size is explained by desiccation tolerance. Larger body size decreased the surface-to-volume ratio, and this is correlated with increased water content 90 . Selection for increased desiccation tolerance showed increased wet weight and lipid content in D. melanogaster 91 . Such proving's support the relationship between size and desiccation tolerance moreover humidity. However, our results are inconsistent with this pattern. The study with field collected Anopheles albimanus populations have exhibited a positive correlation between wing size and relative humidity 92 which is consistent with our findings. The inconsistency among the studies with humidity and body size point out the need for further work to understand the relationship between humidity and also precipitation and body size. It is known that humidity has a large impact of some traits in D. melanogaster [93][94][95] and that rainfall is associated with adaptation in D. melanogaster natural populations 96 . Moreover, our results reveal that humidity and precipitation have a greater impact on males than females where temperature selection acts the same direction in both sexes. This sex-specific effect for precipitation and humidity could be explained by the selection in nature on body size in the context of fitness advantage. Body size is a sexual dimorphic trait but the genetic basis of body size is shared in both sexes. But it is known that these genes have diverge effects by sex on body size in D. melanogaster 97,98 . Like body size, lifespan and ageing share also a same genetic basis in both sexes but sex-specific effects of lifespan were recorded for D. melanogaster 99 and D. simulans 100 . Archer et al. 100 demonstrated that the sexual and natural selection have sex-specific effects in D. simulans where males evolve grater baseline mortality than females under high temperature conditions 100 . Natural selection may be a determinant sex differences for wing size in Yeşilöz population where precipitation and humidity favors larger size in males may due the differences in fitness outcomes.
Populations in seasonal fluctuating environments receive multiple environmental cues and must deal with this heterogeneous environment to survive and reproduce. A growing body of literature shows that this situation is resolved through rapid adaptation in D. melanogaster populations. This study shows that body size, an important fitness component, is one of the changing phenotypic traits across seasons. Our results suggest that the effects of environmental variables that generated variation in body size among populations such as latitudinal clines, act as a selective pressure in a different manner in seasonal variation within population. Long-term monitoring of a population in its natural habitat and quantitative measurement of its responses to seasonal environmental changes are important for understanding the adaptive response of population to temporal selection. Future